%Beta Robust Median estimation
function [beta]  = BetaRobustMed(rStock,rIndex)
    [n,m] = size(rStock);
    X = [rIndex];
    y = rStock;
    %Calculates the Least Median of Squares (LMS) simple
    % p is the number of parameters to be estimated
    p = size(X,2);

    % The "half" of the data points
    h = floor(n/2)+floor((p+1)/2);

    %All the possible combinations of m+1 m-dimensional points
    C = combnk(1:n,p);

    rmin = Inf;
for i = 1:size(C,1)
     for j = 1:p
         A(j,:) = X(C(i,j),:);
         b(j,1) = y(C(i,j));
     end
         if rank(A') == p
            % Calculate  SLOPES
            c = (A'*A)\A'*b;
            est = X*c;
            r = y-est;
            r2 = r.^2;
            r2 = sort(r2);
            rlms = r2(h);
            if rlms < rmin
                rmin = rlms;
                beta = c;
            end
         end
end
